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ABSTRACT 

The  Kalman  Filter  is  a  widely  used  procedure  in  tracking  algorithms.  When 
normality  assumptions  are  violated,  the  Kalman  Filter  performance  tends  to  degrade. 
In  this  thesis  a  new  procedure  is  introduced  for  accommodating  non-normal  properties 
of  measurement  error  distributions.  The  procedure  is  developed  for  the  muiti-observer 
situation.  Simulation  experiment  results  are  presented  and  numerical  comparisons  are 
made  between  the  Kalman  Filter  performance  and  that  of  the  new  procedure. 
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I.  INTRODUCTION 

Tracking  algorithms  are  widely  used  in  modern  technology.  Perhaps  the 
most  commonly  used  algorithm  is  the  Kaiman  Filter  (KF)  (cf.  [Ref.  ij  and  [Ref.  2]  ) 
which  is  based  on  the  following  assumptions: 

a)  The    target    movement    model    is    known,    and    has    normally    distributed 
fluctuations. 

b)  The    sensors    tracking    the    target    provide    unbiased,    normally    distributed 
measurements  of  current  target  ''position". 

If  the  assumptions  are  satisfied.  KF  is  the  best  algorithm  possible.  In  practice, 
however,  the  normality  assumption  may  not  be  justified,  and  the  KF  performance  is 
degraded.  A  number  of  modifications  of  the  KF  have  been  made  by  different  authors 
based  either  on  heuristics,  analytic  derivations,  or  both  in  order  to  improve  KF 
performance  in  different  situations  (cf.    [Ref.  3]  and  [Ref.  4]  ). 

In  this  thesis  a  tracking  procedure  is  developed  for  the  following  one-dimentional 
basic  (and  classical)  model: 

a)  9/I  +  1  =  dn  -  eB+1;  en  +  1~N(0,x2);      «- 0,1,2,. ..« 

b)  yn+lj  =   9*+l  +  *n+lj      *n+U~*dffy       J=l>2<3>"m   «  =  0,L2,...QO 


where 


^njr  i  is  the  target  position  after  time  step  n 


v„_  i  •  is  the  measurement  of  target  "position"  obtained  bv  sensor    at  time  step 

!-   I. 

In  the  model  the  target  oerforms  a  random  walk  with  normally  distributed  steps  with 


'  -3 


ih 


mean  zero  and  known  variance  t  ,  and  the  /  sensor's  measurement  errors  are 
unbiased  and  t-distnbuted  with  known  oarameters  (d.,G.)  where  d.  is  the  decrees  of 
freedom  and  a.  is  the  scaie  parameter  of  the  t-distribution.    There  are  m  sensors. 

This  model  violates  the  classical  KF  assumptions  by  allowing  thick-tailed 
(outlier- prone)  measurement  errors.  This  model  is  aesigned  for  the  situation  where  m  is 
small,  and  it  isn't  possible  to  make  use  of  the  central  limit  theorem  effect. 

The  mathematical  tractability  of  the  model  is  complicated  by  the  fact  that  there 
is  no  simple  analytic  form  for  an  estimate  of  location  of  the  target  with  Student-t 
measurement  errors. 
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In  Chapter  II  two  possible  ways  to  derive  the  classical  KF  equations  are  given. 
These  derivations  will  help  to  motivate  the  procedure  suggested  in  Chapter  III.  Some 
practical  considerations  regarding  the  implementation  of  the  procedure  will  be  given  in 
Chapter  III  as  well. 

In  Chapter  IV  the  simulation  experiment  implementing  the  procedure  derived  in 
Chapter  III  is  described  and  numerical  performance  comparisons  of  the  KF  and  the 
new  procedure  are  presented. 

Chapter  V  contains  final  remarks  about  the  new  filter  and  further  topics  for 
future  development. 

Appendix   A   contains  tabulated   results    of  simulation  experiment. 

Appendix  B  contains  details  on  simulation  experiment  implementation. 

Appendix  C  contains  graphical  presentation  of  the  results  from  Appendix  A. 
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II.  KALMAN  FILTER 

The  basic  model  for  the  classical  Kalman  filter  is  that  the  target  moves  according 
to  a  random  walk  with  normally  distributed  step  sizes.  !  he  measurement  errors  are 
also  normal  and  statistically  independent  from  target  "position'",  i.e.; 

ai       0    .  .  =  0    fg    ,  ,;  c,  ,~N(0,T2);       «  =  0,1,2...  .°o 

/i  -t-  1  n  /r  -r  [ '  n~r  \  v    '      ' '  '    '    ' 

There  are  three  alternative  ways    to    obtain    the    KF  equations: 

1.  the  maximum  likelihood  estimation  approach. 

2.  Tie  minimum  variance  unbiased  estimation  approach. 

3.  the  minimization  of  the  mean  square  of  estimation  crroi  ;0     -  6  ). 
Extension  and  synergism  of  the  first  two  methods  will  lead  us    to  the  procedure 

described  111  Chapter  III. 

A.       MAXIMUM  LIKELIHOOD 

A  A 

Let  8  be  the  estimate  of  8  after  time  step  n  and  C  be  its  variance.  Since  8  is 
normullv1  distributed    with    mean  8    and  variance  C    ,    the  conditional  likelihood  of 


8    1  -  iiiven  8    and  v    ,  ,  , y    ,  .         is: 

/it  I  -  1  -  n  t  1,1  •      •'/i+l  ,m 


a  +  1        i  1      1  +■  1 .;        n  +■  1 


rt-t-l'    Tya-t-l,r      '•/n  +  l,m  '    ' 

7  =  1 


where  i;Q  =  !?,,(//)  =  Cn  +  zz 

and  i'.=  cr.2-  Var(y.)     (/'=  1 m). 


In  the  KF  case  normality  is  preserved  at  each  estimation  step  as  can  be  easily 
verified. 
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Finding  the  maximum  for  this  likelihood  is  equivalent  to  finding  minimum  of  the 
exponent, 
i.e. 

m   (8   ^.-y.)2\ 

min(      '  '      )  =  Q(6    +1)  (2.2) 

(dropping  the  first  subscript  on  the  j't:+['s  and  using  }\]-^n  for  a  compact  notation). 
Expanding  the  squares  we  obtain 


Q(9      ,)  =  82       ''--20         '"  -i+^-i  (2.3) 

j=()    j  J=Q    j        7=0    j 


This  simple  quadratic  expression  is  minimized  by  letting 


A 

8  m     v 

^    —  "  +       .  

a  /  =  0    /  i  y  J 

C/T.l) 


v  1        _L_  f  v  I 


Since  (K.\)  is  a  linear  equation  we  may  immediately  obtain  the  variance  of  the 
estimate. 

/     m    y  *     V 

iK.2)  C  =   Vane       J  =  ' —    =    =    . 

i  + 1  a-t- 1  2  1  m 


Equations  (K.l)  and  (K.2)  are  the  only  ones  needed  to  implement  the  KF. 

b:     minimum  variance 

A 

Given  yn  =  '3ri    and  y-=yn+[        /'=  1 m  ;     it    is   desired    to    get     an    unbiased 

minimum  variance  estimate  (UMVE)  of  0    ,  ,  . 

n  ■+"  1 

Since  yn.;',,...,y  „  are  all  unbiased  estimates  of  0    ...  anv  linear  combination 

H  =  a.v;'n-«,j',  +  ...+  a   y  (2.4) 

constrained  by  Ya.=  1  will  also  be  an  unbiased  estimate  of  0    .  .  . 

J     "— '     J  /IT  1 
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Since 

Vax(H)-aqV«iV»-  +  aA,  •  <2-5) 

finding  the  UMVE  is  equivalent  to  solving  the  following  problem: 

,.         i  •>  •> 

mm  v  =  <l,~vn  +a,"v.  -r  ...  4- a    *v 

U     0         11  m     m  (7  f\\ 

0  1  tfl 

Given  a  Lagrange  multiplier  X,  the  Lagrangian  equation  may  be  expressed  as: 


i£  =  an-vn  +  a,2v,  + ...  +  a  -v    -  X(an  fa.  +  ...  +  or  -  I).  (2.7) 


After  taking    partial    derivatives    and  setting  them  equal    to    zero  the  resulting    system 


2a0v0-X-0 

2a,  v.  -  X,  =  0 


2a   v    -X  =  0 

m  m 

-  a„  -a,  -  ...  -a    =  -  1, 

U  1  <H 


;2.s) 


Multiplying  each  equation  by  the  corresponding   l/(2v.)  and  adding  all  the  equations 

together  results  in: 


V   =    (2.9) 

v  1 


Substitution  then  leads  to 


u 
a-(  -    — —  (2.10) 

J  m 

-   I 
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So,  once  again  the  final  equations  are  equations  (AM)  and  (K.2) 


IK. I)  8         =  H  = 

■1+  I 


V  J. 

-    i 


(K\2)  C      ,    =  mm  V  - 

a-t-  I 


V  1 


Note  that  the  minimum  variance  approach  does  not  require  normality  assumptions, 
but  normality  assumptions  make  it   equivalent  to  the  maximum  likelihood  approach. 
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III.  W-FILTER 

Now  to  return  to  the  basic  model,  target  motion  is  a  random  walk  with  normal 
step    size,   but    the   m   sensors   have   independent   Student-t   distributed   measurement 
errors. 
The  basic  model: 


a)       0    ,  ,  =  0    4-  c    ,  ,; 

n  -r  I  n  «ti' 


c  4.,  ^Ni'O.-r);   «  =  0,1,2... .oo 


b)       v4.l=0.1-rd,l;      o,.    ~t(c/,(7);   /=  l,2,3,...m   «  =  0,l,2,...oo 

It  is  assumed  chat  our  procedure  is  producing  unoiased  approximately  normally 

A 

distributed  estimates.    This  means  that  after  time  step  n,  8    is  approximately  normally 
distributed  with  mean  0     and  variance  C    .     This  is  the  same  assumption  as  used  bv 

a  n 

A 

[Ref.  3]  and  [Ref.  4],    Based  on  this  assumption  we  need  to  construct  9^  +  1  and  C  +  t  . 

A.       THE  BEST  0    .  , 
I.  The  Criterion 

Based  on  the  assumption  above,  the  conditional  likelihood  of  9    ,  ,  given  9 
and>-;i+1  1,."U;rt+l  m  ^as  the  form: 


Lie     |e 


(0 

1     'l"1- 

A 

2 

2 

"     c 

9 

r  t" 

.71 

•n 

7  =  1 

1  + 


ctd ) 


v  -8  i 


d.  +  l 

J 


(3.1) 


where  c(^.)  is  some  constant  depending  on  the  degrees-of-freedom  parameter  d.  of  the 
/  sensor.  Although  it  is  possible  to  look  for  the  maximum  likelihood  point,  taking 
that  approach  has  the  following  two  interrelated  potential  problems,  first  there  may 
be  multiple  local  maxima  present.  In  an  extreme  case  there  may  be  as  many  as  m+  1 
local  maxima.  (Sec  [RcL  A\  for  local  maxima  conditions  analysis  when  m=  1.)  Second, 
in  the  multiple  maxima  situation,  the  global  maxima  may  be  tail  and  skinny,  not 
having  much  likelihood  to  support  it.  One  would  not  like  to  commit  all  belief  on  the 
basis  of  such  evidence. 
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To  deal  with  these  problems  use  the  following  (smoothing)  approach  will  be 

A 

used:  instead  of  choosing  61,  to  maximize  the  probability  of  having  9^  _^_  j  in  an 
infinitesimal  probability  element  dO  around  Qn+y  ,  the  estimator  ^nJrl  will  be  chosen 
so  as  to  maximize  the  probability  of  having  0jt[  in  an  interval 
[6  +  ,  —  w,Q  +  1  +  w],  the  later  interval  having  a  finite  (usually  non-zero)  predetermined 
length  2w.  In  other  words  the  concern  is  not  having  an  estimate  exactly  "on  the  spot" 
but  rather  having  it  within  distance  w  from  true  value  of  0^+  {  . 

In    mathematical   terms    the   approach    will   be    to    find    the   solution   to    the 
equation: 

C  9       +  w 
maxKQ        )  =   max  \  Ui$    y  y         JdZ 

(w.i.a)  \+rw 

n  ■+•  I  n.  +■  1 


In  the  equation  [IV. I. a)  w  is   serving  as  a   tuning   parameter  and    has   the 

following  effect  on  the  solution: 

When  vv  =  0  solving  the  equation  (W.l.a)  is  equivalent  to  finding  maximum 
likelihood  point. 

When  vv>0  and  finite,  equation  (W.l.a)  tends  to  down-weight  skinny  peaks 
more  than  fat  ones. 

When  w  is  large  enough,  equation  (W.l.a)  will  have  a  unique-  solution. 
However  it  is  generally  neither  required  nor  optimal  to  use  very7  large 
values  for  w. 

When  w  approaches  infinity,  any  value  of  9  +  .  will    satisfy    the  conation. 

The    idea    of  using  a  non-zero  w  has  additional  appeal  since,  in  some  practical 

situations,  occurrence  ot  a  large  tracking  error  may  degrade  sensor  performance  for 

consecutive  measurements.    (Such  sensor  dependence  on  tracking  performance  is  not 

reflected  in  the  model.) 

2.  The  Solution  Technique 

Finding  a  solution  lor  equation  (W.l.a)  directly  involves  an  exhaustive  search 

with  numerical  integration;  not  a  very  exciting  prospect  at  first  sight.  However  it  has 

the  advantage  of  algorithmic  simplicity  and  guaranteed  global  maxima. 


"Since  the  likelihood  function  is  positive  and  asymptotically  approaches  zero 
when  G/;+  [  approaches  ±  00  this  statement  can  be  easily  proved. 
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a.  Analytic /iterative  approach 

Alternatively,  we  may  try  to  find  an  analytic  solution  to  the  equation 
(IV. I. a)  by  taking  the  derivative  of  /(Gn+  j)  with  respect  to  0  +  ■  and  setting  it  equal 
to  zero. 


9       ,  +  w  -  8    f 

I      n  +  I  i 


dI(B 


i-t-r 


C      4-1' 

1 


a'9 


i  +  I 


•n 


dd  ) 


j  -i 


[^ 


V  ,      —  'V  — 

"  n  +  1 .; 


cf.  J 


(3.2) 


19  -u;-8    )" 

1     i  +  I  i 


2  „         2 

C    +  t 

i 


•n 

7  =  1 


ci<i. 


d  .  +■  i 


/  y      ,  .  +•  w  —  8    ,  ,  v     i    I     2 

/  -1+  1./  ii+l  L 

1+ I  — 

a.  /  d.\ 


Alter'  setting   the   derivative   equal   to   zero   and     performing  cancellations 


<9       ,  +  uo-Q 
I     n  + 1 


C    +V 

1 


•n 


:    -I 


14- 


'   '  7 

y      .  .  —  w  —  9  i 

"1  +  1,;  i+l 


/     d. 


A  A        7 

(9       ,  -uo— 9   )*" 

li  +  l  i 


(3.3) 


■n 

•n 


i  + 


y       ,    +  w  —  0 

y  i  +  1 .;  n  t-  1    \       1 
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Next,  taking  logs  of  both  sides  of  this  equation  results  in: 
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Multiplying  through  by  (  —  1),  expanding  the  squares,  and  pulling  out  the  common 
denominator  of  the  logarithm  gives: 


92      +  uT  +  02V20      j.v-20  iv-2Q      ,9 
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Collecting  terms  and  rearranging  leads  to 
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And  finally, 
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Equation  (IVA.b)  is  the  first  order  condition  for  a  local  maximum  for 
equation  (IV. La).  But  perhaps  a  more  natural  way    to  check  the  first  order  conditions 

A  A 

at  point  0/i+  j  is  to  see  if  the  likelihood  function  takes  equal  values  at  points  0^+.  .  -  w 

and  0    .  ,  +  w. 

n  +  1 
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When  w  approaches  zero,  equation  (IV. I. b)  approaches  the  form: 

-  m    (d.+  l)(y     .,  .-%    x1) 

9         =9     +  (C  +x2)  V  -J- n  +  lJ  ^  +  \  0.7) 

J      l     j    j       J  n-rlj         n-i-l 

(The  right-most  term  is  the  first  derivative    with    respect    to  9  4.1  of  the  logarithmic 

function.) 

This  is  exactly  the  first-order  condition  for  a  maximum   likelihood  parameter  value. 

Since  9  +,  appears  in  the  logarithmic  term  of  (IVA.b),  equation  (IVA.b) 
is  difficult  to  solve  directly.  One  solution  procedure  is  to  use  {IVA.b)  in  a  recursive 
manner  to  obtain  6,  _t_ ,  .  The  empirical  experience  with  this  approach  indicates  that 
unless  :v  is  relatively  large,  direct  application  of  equation  (IVA.b)  does  not  converge  in 
most  cases. 

b.    Integral  evaluation  approach 

Instead  of  using  a  sophisticated  iterative  procedure,  the  simple  exhaustive 
search  procedure  described  below  will  be  used. 
Define 

A 

0    .    =  mitiiQ  .y    .  ,  ,  ,...,y    ,  ,     ) 
mm  v    ny/n-r\.V     ^n+\.m'  (IS) 

0        —  maxlv  .v    ,  ,  ,,...,}'  _i_ ,     ) 

max  nvn^  1,1'     ^n+\,nr 

A 

The  interval  T=iO    .  ,9       1  includes  all  reasonable  candidate  points    for  9    .  .  .    It  is 

1    mill'    max'  [  «t| 

possible  to  evaluate  equation  [IVA.b)  at  a    finite  number  of  points    (say  20)  from  the 
interval  T  in  accordance  with  the  resolution  required,  and  take  the  one  with  smallest 

A 

violation  of  equation  ( IVA.b)  as  the  estimate  9  j_ ,  . 

t  i~  1 

The  approach  of  picking  the  point  with  the  minimum  violation  of  equation 
(IVA.b)  does  not  guarantee  even  locai  maxima.  However  in  some  of  the  simulations 
that  were  performed  it  proved  satisfactory  unless  there  were  many  cases  when  the 
likelihood  function  iiad  large  "fiat  spots"  and  iv  was  small.  This  generally  happens 
when  the  measurement  variances  are  large  and/or  the  number  oi  observers  m  is  large. 

A  procedure  based  on  evaluation  ot  first  or  higher  order  derivatives 
will  be  computationally  costly  and  in  general  will  be  more  difficult  to  solve  because  oC 
local   maxima   and  "Hat  spots"  in  the  likelihood  function. 


The  "jumps"  between  iterations  may  be  very  large  and  the  estimate  may  reach 
totally  unreasonable  values. 
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Therefore,  the  safest  procedure  for  finding  the  solution  of  equation 
(IVA.a)  is  exhaustive  search  over  the  interval  T  with  numerical  integral  evaluation. 
This  is  the  procedure  used  in  the  simulation  studies  reported  in  the  next  Chapter.  In 
the  simulation,  search  and  integral  evaluation  were  performed  in  a  single  DO-LOOP 
(See  Appendix  3  for  details  on  implementation).  The  simulation  was  nor.  written  with 
computational  efficiency  in  mind,  nevertheless  the  computational  burden  did  not  reach 
an  unacceptable  level. 


B.        VARIANCE  OF  0    , 


*+1 


In  order  to  keep  the  filter  going  once  6,  +  .  has  been  obtained,  it  is  necessary  to 

compute  its  variance  C  +[  .    Since  equation  {IVA.a)  is  implicit  and    equation  (IVA.b) 

is  not  linear,  there  is  no  natural  analogy  to  equation  (A'.2). 

1.    Linearization  approach 

A 

One  possible  approach  to  calculate  the  variance  oi  (i  +.  is  to  approximate 
equation  (IVA.b)  with  a  linear  equation  obtained  by  a  first  order  Taylor  expansion  of 
the  logarithmic  terms  around  8-  ,  ,  =y  .  ,  for  each  /.  Cancelling  and  rearranging 
terms  results  in  the  following  approximation  to  equation  (IVA.b): 


m 

iWl.c)  9      ,   =K~[{  9     +  (C  +t2)y  k.y      ,.) 

n  + 1  n  n  —       ; "  n  +  1  j 

7  =  1 

where 

d  4-1 
*     =    — ;  (3.9) 

i  a ".  -  iv " 
j   j 

m 

K  -  1  +  (C   -t-t2)  y   k  (3.10) 

n  —       j 

; = l 

This  leads  to: 


m 

C        =  K~2{C    +  (C  +  iV  V  k2  v    )  (3.11 

n  +• 1  i  n  —      j     j 

7=1 


where  v.=  Var(y    ,  ,  .). 
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Empirical    evidence  indicates  that  in  practice    this    approach  performs  very 
badly,  giving  values  for  Cn  which  are  far  too  small.  The  Cn  depends  only  on  n  and 
converges  to  a  value  for  large  n. 
2.  Minimization  approach 
a.     The  Idea 

Having    rejected    the    above    linearization    approach,    a    different    one    is 
introduced  by  the  following  motivation. 

The  "best"  estimate,  i3  ,  ,  ,  is  some  combination  )f  9  y  _i_ ,  ,-,...&  _i_  , 
In  the  case  of  the  KF  it  is  possible  to  pick,  a  linear  combination  having  minimal 
variance  and  remain  consistent  with  maximum  likciihood  criteria.  In  the  WF  case 
this  consistency  no  longer  holds,  since  generally  the  estimate  obtained  by  finding  the 
solution  to  the  equation  (IV. I. a)  will  differ  from  the  minimum  variance  estimate.  But  in 
this  case  it  :s  still  possible  find  a  minimum  variance  linear  combination  constrained 
to  vield  the  "best"  faireadv  known)  value  tor  0    .  ,  . 

This  means  that  in  order  to  obtain  an  estimate  ol  variance  C    .  ,  we  have 
to  solve  the  following  constrained   minimization  problem: 

The  svmbolic  formulation  would  be 


nun  V  =  (X,, "v., -r  (Jt"y,  +-...  — U    ~v 
u     U         1      i  '7i     m 

S.T.    a,,0   +a,i-       ,  .4-...  +  a   v    .  .      =0    .  , 

fL  +  a,  4-  ...j- a    =1 


,3.12) 


lo  '  "i  ' 


m 


where 


on  =  C  -c2  (3.13) 

0  a 


o     =    Variy       ,  .  )     = J. 


;=l.2....,m  (3.14) 


solution. 


For  mult:  derivation 

If  m=  1.  Lhen  no  minimization  is  necessary  and  the  problem  has    a  unique 

A 

For  notationai  simplicity  the  first  subscript  will  be  dropped,  and  let  j,'0  =  8 
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A       A 

and  9  =  0^+1  .    Let  Xj  and  X.,    be     Lagrange  multipliers.  The  Lagrangian  formulation 
is  then: 


££=  V  aV    -  \t(  Y  a^-9  j  -  AJ^   >.  a.- 


7  =  0  7  =  0  ;=0 


Taking  partial  derivatives  with  respect  to  a.  and  X.  ,  X.,  to  obtain    a  system  of  m+3 


linear  equations  with  m+  3  unknowns  gives: 


2v0aQ  -  \{yQ  -  X2  =  0 
2i'ja,  -  "k,y.  -  X->  =  0 


2v  a     -  X,-.'     -  X0  =  0  (3.16) 


7=0 


;=0 


Multiplying  each  one  of  the  first  m  +  1  equations  by  the  corresponding  1/v.   and  adding 
all  together  results  in: 

■n     v  •*      , 

_v    V  J._K    V  1  =   _2  (3.17) 

i  —  *>  

;=0    7  7  =  0    ; 


Since  ^o.=  1. 

Multiplying  each  one  of  the  first  w+  1  equations  bv  corresponding  y.,'v.  ,  the  (m  +  2) 


equation  by  (  —  2)  and  adding  all  together  produces: 


VI      y  Tt      v 


_\  y  -i_x  y  -i  =  -29  o.is) 

7=0    7  7=0    J 
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After  changing  signs  a  system  of  two  linear  equations  follows: 


XlB+X2C=2 


where 


A 


Applying  Cramer  s  rule  gives 
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Substitution  vields 
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(3.19) 


28C  -  25 

V   = (3.20) 

AC  -  B2 


2A  -  28/3 
\.„  =    .'3.21) 


a    =  ■ l/'  =  0,l,2,...,m)  (3.22) 


(3. 


v         n  rr,       v 

^iiy  LJ  -jl\ 


Note  that  the  denominator  is  equal  to  zero  in  case>,0=>,]  =..:=j-     .  which  is  a  highly 
improbable  event. 
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Since  the  Hessian  matrix 
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is  positive  definite  (  v.>0]  for  ail  y),  the  solution  is  a  global  minimum. 

J 

The  filial  equation  lor  the  variance  estimate  is: 


,     m     v"    ti  /     m     v       ~ )  m 

JyiLy.L.Jyy2\\    v 
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4. 
5. 


Now  the  equations  (WA.a)  and  ((T.2)  contain  all  that  is  necessary  to  implement  the 
WF  procedure. 

C.  THE  WF  PROCEDURE 

1.  Find  0    •    =//ii//iO  .y    ,  ,  ,,....y    .   ,      ),      9        =  ma.v(9  .y    ,   ,  ,,...,>'    ■  ,     ). 

nun  fr/i+l,!1     u/iriy  max  v    /zv/r-r  1,1'     ^jit-  l,/n' 

2.  Fix  the  actual  search  resolution  R  (depending  on  T==|G        ,9       !). 

v       f  '    min        max" 

3.  Set  up  the  grid  with  resolution  R  over  the  interval  (0    .    —  w,  0        4-w). 

1  v    mm  '      max 

(Let  k=-u>/R) 

Evaluate  the  likelihood  function  at  each  grid  point. 

Find  the  9>;+  [  by  picking  the  "best"'  grid  point  solution  ot    equation  (  WA.a)  in 

interval  T  using  sum  of  2k +  1    adjacent    likelihood  ('unction  values  lor  integral 

approximation. 

-Sue  Appendix  B    for  details  on  implementation). 

6.      Compute  C    .  ,  bv  equation  t  W.2). 

D.  PRACTICAL  CONSIDERATIONS 

A  number  of  issues  related  to  implementing  the  WF  procedure  on  a  digital 
computer  arc  listed  below.  These  considerations  were  used  in  in  simulation  described  in 
the  next  Chapter. 

1.   Grid  setup 

The  exhaustive    search    over    the    interval  T  =  iO    .  ,0       ]  mav  be  performed 

L    mm-    max1  -  r 

over  a  pre-specified  number  of  points  to  be  checked  or  using  a  pre-specified  search 
resolution. 
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The  first  approach  is  wasteful  when  T  is  small  and  the  second  approach  is 
wasteful  when  T  is  large.  Thus,  it  would  seem  reasonable  to  specify  both  a  desired 
resolution  and  an  upper  bound  on  the  maximum  number  of  points  to  be  checked, 
which  remain  constant  over  time. 

2.  Very  dense  observations 

In  the  case  when  T  is  very  short  [here  mav  be  two  complications: 

i)  Because  or  resolution  problems  9  _,_  L  may  take  exact  value  of  9  .  or  Q;riax  ■ 
In  this  case  equation  (IV.2)  may  force  C  +,  equal  to  the  corresponding  v., 
which  is  not  right. 

Z)  Computing  C  .  ,  by  (IV. 2)  with  y.  values  verv  close  together  mav  yield 
underflow    ana      division      by  zero. 

In  order  to  prevent  these  problems  the  following    modification  is  suggested 

and  was  used  in  the  simulations:  if  T  is  very  short  (say   less   than   3    times  resolution) 

A 

set  9  _l.  L  to  be  equai  to  midpoint  of  T   and   compute  C  +  {  using  (K.2). 

3.  "Very  small"  numbers 

When  m  is  iarge,  the  likelihood  function  takes  very  low  values.  To  avoid 
underflow,  simply  multiply  it  by  a  iarge  constant.  The  constant  that  was  used  in  the 
sunulation  is  1024'". 

One  more  modification  was  implemented  to  avoid  arithmetic  with  very  small 
numbers.  If  the  exponent  of  the  likelihood  function  first  term  was  less  than  -  100. 
the  likelihood  function  was  set  to  be  equai  zero.  The  meaning  oi  the  modification 
is  that  the  procedure  considers  any  target  step  larger  than 
^7  x  (standard  deviation)  to  be  a  practically  an  impossible  step.  (One  must  be  careful 
with  such  modifications  if  large  target  jumps  are  possible,  as  could  occur  if  the  model 
involved  Student-t  9-iumps). 
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IV.  SIMULATION  RESULTS 

A.       SIMULATION  DESCRIPTION 

A  simulation  experiment  was  performed  on  an  IBM  3033  and  4381  at  the  Naval 
Postgraduate  School.  The  programming  language  used  was  FORTRAN-77,  and 
single  precision  was  used.  Normal  and  Chi-square  random  numbers  were  generated 
with  the  aid  of  IMSL  library  routines.  Student-t  random  numbers  were  generated 
from  normal  and  Chi-square  random  variables  with  appropriate  degrees  ol  freedom 
(see  [Ref.  5]  ).  Tracking  error  statistics  were  computed  using  the  HISTGP  subroutme 
from  the  NONIMSL  library  at  NPS.  Figures  were  generated  using  an  IBM 
experimental  graphics  package.  GRAFSTAT. 

For  each  simulation  replication  the  tracking  sequence  was  performed  for  100 
steps  ,  (i.e.  n  +  1  =  1,2 100).  The  number  of  replications  was  1000  in  ail  cases. 

The  KF  and  WF  procedures  (with  different  values  for  w)  were  carried  out  using 
the  same  random  numbers.  (With  normally  distributed  target  steps  and  Student-t 
distributed  measurement  errors). 

A 

Statistics       were       collected       on       tracking       errors      9    .  ,— 0    .  ,       for 

°  n  "t-  1         n  T 1 

«+•  1  =  1,25,50,75,100. 

The  simulation  was  performed  for  6  cases: 
Case  1  -one  observer,   cr,  =  1 
Case  2  -two  observers,    c\  =  <3\  =  1 
Case  3  -three  observers,    ff,  =  d.,  =  <5\  =  1 

1  2  j 

Case  4  -three  observers.   <r;  =rcri=  1:      (7,  =  3 

L  2  J 

Case  5  -five  observers,   a,  =  d0  =  d,  =  d4  =  0"5  =  1 
Case  6 -five  observers,   d,  =<t^  =  <t,  =  1;      d.  =  (y-=3 

1  2         j         '  4         5 

In  all  cases  the  measurement  errors  of  the  observers  have  the  Student-t  distribution 
with  d=  3  degrees  oi  freedom  and  scale  <y;,  and  the  target  has  normally  distributed 
step  sizes,  with  standard  deviation  T=  1.  In  all  cases  the  WF  procedure  was  performed 
for  values  vv=  0,1.2,3.  Excluding  case  1,  the  desired  grid  resolution.  DR,  was  chosen  :o 
be  0.1.  In  case  1,  the  desired  resolution  was  equal  to  0.05.  Case  1  is  most  susceptible 
to  problems  resulting  from  a  very  short  interval  T  described  in  the  previous  chapter. 
In  all  cases  the  number  of  grid  points  on  the  interval  T  was  limited  by  upper  bound 
of  50.  (See  Appendix  B   for  details  on  resolution  and  grid  setup.) 


27 


After  fixing  a  grid,  the  likelihood  function  was  evaluated  at  each  grid  point  over 
the  interval  of  length  T  +  2w.  The  integral  was  approximated  by  the  sum  of  2k  +1 
values   of  the  likelihood  function  with  k  being  the  largest  integer  not  exceeding  w/R. 

B.  WF  PERFORMANCE 

Figures  C.l  through  C.6  display  comparisons  of  the  standard  deviations4  ot 

A 

the  tracking  error.  9._,  ~~ ®n-i-[  •  for  the  KF  and  the  WF  for  cases  1  through  6, 
respectively.  In  all  the  cases,  the  WF  performance  is  considerably  better  than  the  KF 
performance.    See   Appendix  A  for  tabulated  results. 

Figures  C.7  and  C.S  display  WF  tracking  error  standard  deviations  with  different 
values  of  w.  The  optimal  w  is  somewhere  in  the  vicinity  of  1  or  2.  but  sensitivity  to  the 
exact  value  of  w  is  low.   See  Appendix  A  for  exact  numbers. 

C.  SENSITIVITY  ANALYSIS 

For  sensitivity  analysis  purposes,  the  tracking    sequence   was  performed  on  each 

case  with  the  same  steps  and  measurements,  but  the  WF  procedure  used  the  following 

(wrong)  assumptions: 

<J.  =  0.5    or    <r;=  1.5    instead    of    ff.=  1.0    or    <T.=  3.0    (underestimating    the 
measurement  variance  by  a  factor  of  4) 

d;=6    instead  of  d.=  3  (underestimating  the  measurement  variance  by  a  factor 
of  2)  J 

r  =  0.5  instead  of  t=  1  (underestimating  the  target  step  size  variance  by  a  factor 
of  4) 

Additional  runs  were  performed  when  target  step  size  was  a  Student-t 
distributed   variable  having  3    degrees  of  freedom   and  variance  equal  to  1. 

In  all  sensitivity  runs  the  WF  used  w=  2  and  the  desired  resolution  was  equal  to 
0.1.  Also,  the  maximal  number  of  points   in  the  interval  T  was  50. 

Figures  C.9  through  C.14  display  WF  performance  under  sensitivity  analysis 
conditions  for  cases  1  through  6,  respectively.  Clearly  the  WF  procedure  is  very  robust 


All  the  resuits  presented  in  terms  of  sample  standard  deviations,  i.e.  a  square 
roots  of  the  unbiased  estimates  of  population  variances.  Means  are  not  presented 
since,  in  all  cases  for  both  filters,  mean  error  values  were  very  small.  RMS  {square 
Roots  of  Mean  Square  error)  values  were  also  computed  in  the  simulation.  In  all  runs 
inspected  (the  vast  majority  from  ail  runs)  RMS  values  were  were  slightly  smaller  than 
the  standard  deviations  values  (due  to  division  by  N  in  RMS  computation  and  division 
by  N  —  1  in  standard  deviations  computation). 
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with  respect  to  measurement  error  parameter  estimation  and  less  robust  with  respect 
to  major  alternations  in  target  movement  model.  See  Appendix  A  for  tabulated 
results. 

The  strange  behavior  of  the  standard  deviation  plot  with  t-distributed  target 
step  sizes  is  explained  by  the  following  phenomena:  in  very  rare  cases  when  the 
estimate  9  and  actual  target  position  9  + .,  are  extremely  distant  (the  target  makes  a 
huge  step),  WF  starts  to  ''ignore"'  the  incoming  measurements,  regarding  them  as 
outliers  and  sticks  to  the  old  estimate.  In  such  situation,  Cn  is  increasing  steadily  but 
slowly  (by  It  each  time  step  -  due  to  the  wrong  target  model).  After  some  time 
(which  may  be  as  many  as  30  time  steps),  C  becomes  very  large  and  WF  stops 
disregarding  observations  and  the  estimate  shifts  toward  the  actual  target  position. 
(The  last  modification  described  in  the  previous  Chapter  was  harmful  in  this 
situation.) 

Such  undesired  behavior  may  be  treated  by  ad  hoc  adjustments  to  the  WF 
procedure.  For  example,  after  detecting  a  steady  increase  in  C  we  may  arbitrary  set 
Cn  equal  to  a  large  number  and  recalculate  the  last  estimate.  This  kind  of  adjustment 
may  be  proper  when  actually  implementing  WF  if  there  is  concern  about  the  normality 
of  target  step  size.  None  of  these  adjustments  were  implemented  in  the  simulation 
experiment.  It  is  anticipated  that  if  an  outlier-productive  model  of  the  target  motion 
were  used,  the  effect  noted  would  be  automatically  reduced,  but  to  date  no 
investigation  has  been  made. 

Another  comparison  was  made  between  the  WF  and  KF  procedures  by 
performing  a  simulation  experiment  with  normally  distributed  measurement  errors. 
The  simulation  was  performed  for  cases  la-6a.  Cases  la-6a  correspond  to  cases  1-6. 
respectively,  with  normally  distributed  measurement  errors  having  the  same 
variance  as  t-distributed  errors  in  the  original  cases.  Once  again  the  WF  used  w—  2 
and  the  desired  resolution  was  set  zo  0.1  with  a  maximum  number  of  points  in  the 
interval  T  of  50. 

Figures  C.15  and  C.16  display  the  comparison  of  WF  and  KF  procedures 
for  cases  ia-6a.  It  is  evident  that  the  advantage  of  KF  over  WF  in  cases  la-6a  is  less 
than  the  superiority  of  WF  over  KF  in  cases  1-6.  See  Appendix  A  for  tabulated  results. 

In  order  to  see  how  much  measurement  error  assumptions  may  be  violated, 
several  more  runs  were  made  using  the  Cauchy  distribution  ( Student- 1  with  I  d.f.) 
for   generating    measurement   errors    (WF   assumed    3    d.f.     and   the   correct   second 
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parameter  (J).    In  general  the  tracking    was  steady.  However  in  very  extreme  cases, 

A 

when  the  interval  T  was  extremely  long  and    the  old  estimate  0    was  not  one  of  it's 

endpoints   (actual  resolution   in   such   cases  became   extremely   large),  WF   produced 

an  unreasonable  estimate  and  did  not  recover  from  the  huge  error. 

To    overcome     the   problem,   two    alterations   were     made     to     the    simulation 

implementation: 

The  Maximum  number  of  points  to  check  was  set  to  500. 

In  the  case  when  ail  grid  points  on  the  interval  T  corresponded  to  the  integral 

having  a  value  equal  to  0,  the  new  estimate  was   set    to    be  equal  to  o    ,  the  old 

estimate.    The  original   implementation    produced  0    .    for  an    estimate  in  such 

min 

cases.  (See  Appendix  3  for   a    better  approach  by  altering  grid  setup.) 
After  making  these  two  alterations,  the  runs   were   repeated.   The  tracking   was    steady 
with  no  special  problems.    The    results    are  tabulated  in  Appendix  A  and  may  be 
considered  very  good. 

A  selected  number  of  cases  was  chosen  to  perform  a  limited  variability 
demonstration.  Simulation  runs  were  made  with  5  different  pairs  qi  seeds  ( one  for 
target  steps  and  one  for  measurement  errors)  for  1,2,3  and  5  identical  observers  and 
WF  using  w=2.  The  results  are  displayed  in  Figure  C.17.  Generally,  variability 
between  different  pairs  of  seeds  does  not  exceed  the  variability  between  separated 
time  steos  for  a  sinsie  seed. 
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V.  FINAL  DISCUSSION 

The  WF  procedure  has  a  very  general  nature.  It  may  be  applied  directly  to 
any  type  of  measurement  error  distribution.  The  distribution  of  measurement  error 
must  have  finite  variance  and  known  density  (even  in  tabular  form).  The  WF  filter 
performance  with  exotic  distributions,  using  right  or  wrong  distributional 
assumptions,  needs  to  be  investigated. 

(To  apply  the  WF  procedure  to  models  with  non-normai  target  step  size  there 
must  be  an    efficient    way    to    evaluate    the    conditional  likelihood  function,  or  else 

A 

additional  assumptions  have  to  be  made  regarding  distribution  of  the  estimate  9 
before  and  after  movement  of  the  target.) 

When  applied  to  models  with  normally  distributed  measurement  errors 
(and  normally  distributed  target  step  sizes),  the  WF  procedure  (with  any  w)  is 
theoretically  identical  to  the  KF  procedure.  In  practice,  however,  there  might  be  some 
difference  in  performance  due  to  finite  resolution. 

WF  has  a  built-in  tuning  parameter:  the  value  of  w.  The  optimal  value  for  w  is 
definitely  not  zero,  but  the  sensitivity  of  the  procedure  to  the  exact  value  of  w  is  low. 
The  best  values  of  w  seem  to'  be  in  the  range  (t,  2t),  where  "  is  the  standard  deviation 
of  the  target  step  size. 

Instead  of  using  w  as  a  static  tuning  parameter  where  it  is  kept  constant,  w  may 
be  computed  dynamically,  based  on  the  latest  observations.  For  example,  w=T/4 
constramed  to  w>>0.3  and  w<4  and  the  usual  procedure  for  grid  setup  (see  Appendix 
B)  produced  very  good  results  for  all  cases,  but  was  not  superior  to  the  results 
obtained  with  constant  w.  (Dynamic  w  results  were  not  represented  in  the  previous 
Chapter). 

Another  modification  to  be  explored,  if  it  can  be  justified  by  operational  reasons, 
is  use  of  a  non-symmetric  interval  in  equation  (J-Kl.a),  i.e.  performing  integration 
with  lower  limit  9  + ,  —  w,  and  upper  limit   9  + ,  +  w~  . 

Equation  (IV. 2)  may  be  used  not  only  as  part  of  WF    procedure,  but  as  part  of 

A 

any  "good"  procedure  specifying  9^+1  and  looking  for  C  +  1  . 


optimal  rule  for  dynamic  w  calculation  may  be  discovered  by  further  research. 
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The  most  noticable  difference  between  the  KF  and  the  WF  with  respect  to  Cn 
estimation  is  that  in  the  KF  case  Cn  depends  only  on  n,  converging  to  some  value  as  n 
increases,  and  in  the  WF  case  Cn  is  more  dynamic,  depending  on  the  actual 
measurement  values.  The  dynamic  behavior  of  Cn  reflects  the  variable  amounts  of 
information  received  at  each  measurement  step.  This  phenomenon  gives  an  "'inside"' 
performance  measure  on  how  weil  the  filter  is  doing.  None  of  the  cases  checked  by 
simulation  required  corrective  actions  due  to  abnormal  C  behavior  (i.e.  steady 
increase   or  large  jumps),  except  when  target  step  sizes  were  Student-t  distributed. 

A  limited  insDection  of  C  vaiues  calculated  bv  the  WF  procedure  indicates  that 
the  C  vaiues  computed  are  somewhat  larger  than  the  actual  variances  of  tracking 
errors.  It  may  conjectured  that  there  is  room  for  improvement  in  the  WF 
procedure.  On  the  other  hand,  forcing  C  close  to  its  true  value  (known  from  previous 
simulations)  produces  good  but  not  superior  results  to  those  obtained  in  the  regular 
way,  using  equation  (W.2).  This  suggests  that  the  WF  procedure  is  fairly  robust 
with  respect  to  modest  changes  in  the  C  computation  process.  (A  constant  value  for 
C    mav  be  used  if  there  is  sufficient  knowiedse  of  the  environment,  i.e.  the  C    values 

n  -  n 

that  should  be  obtained  are  approximately  known  and  no  large  jumps  in  C   values  are 
expected). 

Extending  the  WF  procedure  to  the  multidimensional  case  does  not  seem  to 
pose  conceptual  problems.  Equation  (W.l.a)  may  be  replaced  by  an  integral  over  a  k- 
dimensional  rectangle.  The  analog  to  aquation  (W.2)  may  be  obtained  by  minimizing 
the  determinant  of  the  covariance  matrix.  However,  the  computational  burden  of 
solving  equation  (W.l.a)  in  the  k-dimensional  case  and  minimizing  the  determinant 
may  be  excessive.  To  actually  implement  the  WF  procedure  in  the 
multidimensional  environment  would  require  an  efficient  way  to  solve  equation  (  W.l.a) 
and   to  minimize  the  determinant  of  the  covariance  matrix. 
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APPENDIX  A 
TABULATED  RESULTS 

This  Appendix    contains    tables     of  tracking    error    standard   deviations.    The 

cables  are  organized  by  observers  configuration. 


■) 


Tabie  1 


o 


One  ooserver. 


Table  2:    Two  identical  observers. 

3.  Table  3:   Three  identical  observers. 

4.  Table  4:   Three  observers,  one  of  whom  is  3  times  less  accurate  than  the  others. 

5.  Table  5:    Five  identical  observers. 

6.  Tabie  6:    Five  observers,  two  of  whom  is  3  times  less  accurate  than  the  others. 
The  following  notation  is  used  for  specifying  target  step  size  and  measurement 

error  distributions: 

•  N(v)  -  normal  with  mean  0  and  variance  v 

•  t(d,s)  -  Student-t  with  d  degrees  of  freedom  and  scale  parameter  ff  =  s 

When    two  types  of  observers  are  involved,  the  parameters  are  given  for  the  more 

accurate  ones. 

All   cases    above   the    double    line    use   precisely   the   same     target    steps    and 

measurement  errors;  only  the  filters  are  different. 

Filter  notations: 

KF  Kalman  filter  using  correct  variance 

WF  w=i  W-filter  using  3  d.f.,  u>=i  and  correct  variance 

WFvv=i,l/2<y  W- filter  using  3  d.f,  w=i  and  assuming  the  scaie  of  the  t  is  :•" 
rather  than  the  true  parameter  <7  for  each  observer 

Wfw-=i,2d  W- filter  using    w=i  and  assuming  that    the  Student-t  degrees  of 

freedom  are  2d  rather  than  the  true  value   of  d 

WFw=  l,ViX  W-fiker  the  correct  measurement  parameters  and  assuming  the 
standard  deviation  of  the  target  step  size  is  V*x  rather  than  the 
true  vaiue  of  t 
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TABLE  1 
ONE  OBSERVER 


Mvmnt 

step 

Measur. 
Error 

FILTER 

1 

Track.  Lensth 

25            50            75 

100 

N(l) 

c(3,l) 

KF 
WF  w=0 

WF  n/=  1 

wf  w= : 

WF  iv  =3 

0.830 

0.842 
0.832 

1.835 
0.837 

1.108 
1.160 
1.096 
1.024 
1.033 

1.079 
1.082 
1.037 
0.998 
1.021 

1.121 
1.176 
1.124 
1.040 
1.060 

1.117 
1.094 
1.048 
0.993 
1.031 

-  '- 

-". 

WFw=2,i-i(T 

WFvv=2.2^ 

WFw  =  2.!:T 

0.846 

0.821 
0.915 

1.017 
1.013 

1.146 

1.003 

0.993 
1.192 

1.042 

1.038 
1.209 

0.993 
0.998 
1.221 

U3..5S) 
N(D 

t(l~l) 

WF  :v=2 

wf  w=: 

0.853 
0.884 

1.632 
1.349 

1.702 

1.445 

1.141 
1.362 

1.125 
1.406 

N(l) 

N£3) 

KF 
»V  t  w  =  _ 

0.357 
0.871 

1.203 
1.250 

1. 166 

1.209 

1.124 
1.185 

1.131 
1.192   1 

TABLE  2 
TWO  OBSERVERS 


Mvmnt 
step 

Measur. 
Error 

FILTER 

1 

25 

"rack  Length 
50           "5 

100 

NCI)       i 

t(3,l) 

KF 
WF  w=0 
WF  w>=  1 

WF  :v=2 

WF  w=3 

0.735 

0.672 

0.670 

1.682 

0.713 

0.891 
0.801 
0.792 

0.787 
0.820 

0.859 
0.757 

0.^46 
0.743 
0.779 

0.888 
o."84 
0.7"""' 
0.777 
0.816 

0.879 

0.-71 

0.766 

'.-69 

0.808 

■?. 

,". 

WFw=2,y2o" 

\VFw=22d 

WFw=21/it 

0.684 

0.683 
0.835 

0.786 
0.796 

0.909 

0.7J6 
0.749 
0.902 

0.776 

0.734 

1.933 

0.767 
0.775 
0.9.17 

t(3..58) 

N(iy 

t(Y,l) 

WF  w—l 
WF  w-2 

0.749 
).813 

'  "°5 
L004 

'  4Q9 
L080 

0.838 
1.006 

0.857 

1.023 

N[l) 

N(3) 

KF 
WF  .•.•=2 

0.750 
1.778 

0  °~>4 
i!  '59 

0.955 
).989 

0.908 
1.946 

0.914 
1.957  ! 
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TABLE  3 
THREE  IDENTICAL  OBSERVERS 


Mvmnt 
step 

Measur. 

Error 

FILTER 

1 

Track  Lensth 

25           50           75 

100 

N(l) 

tQ,l) 

KF 
WF  iv  =0 

WF  w=  I 
WF  \v=  2 
WF  w=  3 

0.663 

0.587 
0.586 
0.600 
0.040 

0.774 
0.682 
0.673 
0.669 
0.694 

0.748 
0.660 
0.653 
0.654 

0.687 

0.746 
0.633 
0.632 
0.645 
0.680 

0.778 
0.667 
0.658 
0.662 

0.707 

■    m 

~y_ 

WFvv  =  2.::  (7 

WFw=2  2d 

WFw=2!J/2T 

0.602 
0.601 

0.778 

O.60S 
0.674 

0.775 

0.65S 
0.659 

0.802 

0.646 
0.652 
0.793 

0.662 
0.673 
0.820 

t(3,.58) 
N(l) 

t(U) 

WF  w=2 
WF  w=2 

0.661 
0.750 

1.073 
0.843 

1.411 
0.890 

0.720 
0.859 

0.709 
0.836 

N(l) 

N(3) 

KF 
WF  w=2 

0.723 
0.745 

0.815 
0.849 

0.7S9 
0.809 

0.794 
0.826 

0.768 
0.796 

TABLE  4 
THREE  OBSERVERS,  ONE  3  TIMES  WORSE 


Mvmnt           Measur.        FILTER 
step               Error 

Track.  Lensth 
1            25           50           75            100 

N£l) 

t(3,l)            KF 

1  WF  w  =  0 
WFw=l 
WFw=2 

WFw=3 

0.721 
0.6^2 
0.641 
0.659 
0.096 

0.887 
0.793 
0.778 
0.768 
0.795 

0.866 
0.766 
0.758 
0.736 
0.782 

0.871 
0.762 
0.754 
0.716 

0.784 

0.874 
0.769 

0.760 
0.729 
0.802 

?:    |     v. 

WFw=2'/2<T 
WFw=2  2^ 
WFvv=2  \ViX 

O.606 
0.659 

0.825 

0.796 

0.775 
0.878 

0.774 
0.756 
0.905 

0.778 
0.762 
0.897 

0.773 

0.768 
0.923 

t(3..5S)                -"-         |WFw=2           0.740       1.079   j    1.455 
SiD               t(l,l)       WFw=2           0.807  |    0.963       0.989 

0.809 
0.937 

0.805 
0.958 

N(l)               NTC3) 

KF 
WF  w=2 

0.770 
0.798 

0.914  I    0.908       0.900   1    0.878 
0.958       0.937       0.948   j    0.914 

1 
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TABLE  5 
FIVE  IDENTICAL  OBSERVERS 


Mvmnt 

step 

Measur. 

Error 

FILTER 

1 

Track  Length 

25           50           75 

100 

N>!} 

tQ,l) 

"IT 

WF  w=0 

WF  n/=  i 

WF  ',<,•=  2 
WF  vy=3 

0.565 

0.476 
0.480 

0.496 
0.547 

0.656 
0.536 
0.535 
0.53 

0.577 

0.616 
0.528 
0.524 
0.530 
0.567 

0.632 

0.534 

0.532 
0.545 
0.588 

0.603 
0.512 
0.507 
0.513 
0.588 

-"- 

WFw=2,1/i<y 
\VF\v=2,2d 

WFa-=2.!:r 

0.508 
0.497 
0.681 

0.555 
0.542 
0.629 

0.546 
0.534 
0.640 

0.560 
0.552 
0.677 

0.535 
0.517 
0.631 

i 

i 

t(  3..5S) 

N(iy 

t(l,l) 

WF  w=2 
WF  w=2 

0.524 
0.671 

1.520 
0.668 

1.396 
0.689 

0.552 
0.683 

0.527 
0.655 

N(l) 

N<3) 

FCF 

WF  w=  2 

0.618 

0.635 

0.657 
0.695 

0.663 
0.686 

0.650 
0.681 

0.655 

0.675 

TABLE  6 
FIVE  OBSERVERS,  TWO  3  TIMES  WORSE 


Mvmnt 
step 

Measur. 
Error 

FILTER 

1 

Track  Length 
25           50           75 

100 

N(l) 

t(3,l) 

KF 
WF  w=Q 

WF  w=  1 
WF  w  =  2 
WF  n/=3 

0.645 
0.57^ 
0.567 
0.586 
0.617 

0.787 
0.647 
0.641 

0.645 
0.671 

0.739 
0.641 
0.633 
0.639 
0.670 

0.751 
0.647 
0.645 
0.654 
0.697 

0.700 
0.599 
0.601 
0.611 
0.649 



• 

\VFw=2.Y*(S 

\VFw=2'.2a 
WFw=2;iiT 

0.606 
0.584 
0.760 

0.669 
0.644 

0.745 

0.657 

0.644 
0.761 

0.677 
0.664 
0.793 

0.636 

0.616 
0.749 

t(3..58) 

M  1 1 

ti  l".  1 ) 

WF  w=2 

WF  w=  2 

0.635 
0.S02 

0.T20 
0.S01 

1.435 
0.830 

0.690 
0.828 

0.630   i 
0.852 

Nix) 

N(3) 

KF 

WF  w=2 

0.694 
0.724 

0.771 

0.824 

0.790 
0.817 

0.782 
0.827 

0.772  ; 
0.803 
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APPENDIX  B 
DETAILS  ON  IMPLEMENTATION 

In  this  Appendix  the  implementation  of  the  simulation  experiment  is  discussed. 
(There  are  many  different  ways  to  implement  the  WF  procedure.)  In  particular  the 
implementation  of  the  grid  setup  and  search  for  maximal  integral  over  interval  with 
length  2w  will  be  explained. 

1.        GRID  SETUP 

Two  parameters  (excluding  w)  control  the  grid  setup: 

1)  DR  -         the  desired  resolution  (usually  DR  =  0.1) 

2)  NMAX  -  the  maximal  number  of  points  to  search  (usually  NMAX  =  50) 

The  first  step  is  to  rind  interval  T  =  [9    .„ ,t)     J,  where 

1  L    mm',  max-1' 

A 

9    .    =mw(0    ',•  _i_ ,  , y    i  -,     ), 

mm  n^n-r  1,1        ^  n~rl,mn 


9        =  maxiQ  .y    .  ,  .......y  _i_  ,     ). 

max  v    «^nTl,I'      vn>i,mJ 

After  finding  T,  the  desired  number  of  search   points,  n  is  computed:    n=(T/DR)+  1. 

If  n>  NMAX,  then  n  is  set  to  be  equal  to  NMAX. 

In  the  next  step,  the  actual  resolution  R  is  computed:  R  =  T/n  . 

After  fixing  R,  the  number  of  points  on  interval  w,  k,  is  computed:    k  =  [w'R].  (If  w  is 

computed  dynamically  its  value  may  change). 

Now  the  grid  is  fixed  to  have  total  number  of  n+2k+l  equally    (R)  spaced  points 

having  n.+  1  points  on  interval  T  itseif  (one  on   each  end). 

The    likelihood  function  values  are  evaluated  for  each  grid    point  and  stored  in  an 

array.  (In  the  same  order  that  the  grid  points   appear  on  the  real  number  line). 

A  more  economical  way  to  set  up  the  grid  would  be  as  follows: 
First,  evaluate  the  likelihood  function  at  the  endpoints  of  interval  T.  If  the 
likelihood  is  positive  at  both  enas.  proceed  as  described  above,  if  the  likelihood  is 
essentially  equal  0  on  one  or  two  endpoints,  then  eliminate  the  corresponding 
observation,  recompute  interval  T,  and  repeat  the  procedure.  (This  was  not 
implemented  in  the  simulation  experiment.) 
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2.       SEARCH  FOR  MAXIMAL  INTEGRAL 

Once  all  likelihood  values  are  stored  in  the  array,  the  search  for  the  maximal 
integral  was  implemented  in  the  following  way. 

Initiation  step:  compute  sum  of  first  2k +1  values  from  the  array.  This  sum  is 
the  intesrai  approximation  corresponding  to  9    .    .    This  is   the  first  candidate  for 

n<    1 

Search  loop:  advance  on  the  array  each  time  adding  a  new  element  to  the  sum 
and  subtracting  the  oidest  one.  if  the  new  sum  is  larger  than  the  largest  among  the 
oid  ones,  take  the  corresponding  grid  point  as  the  best  candidate  so  far  and  record  the 
largest  sum. 

After  n+  1  repetitions  the  procedure  gives  the  best  candidate  among  the  n  —  1  grid 
points  on  interval  T  (and  the  corresponding  integral  value). 

This  implementation  approach  was  chosen  mainiy  because  of  its  algorithmic 
simplicity.  Adding  and  subtracting  REAL  numbers  repetitively  might  introduce 
some  round-off  error.    This  is  not   a  real  problem  in  present  case  because: 

•  50  additions  and  subtractions  -is  not  a  very  large  number. 

•  The  interest  is  in  the  grid  point,  not  the  corresponding  integral  value. 

•  In  close  cases  when  round-off  error  is  reversing,  the  test  can't  be  sure  anyway 
because  of  finite  resolution. 
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APPENDIX  C 
COLLECTION  OF  FIGURES 

In  this  Appendix  seventeen  Figures  presenting  the  WF  performance  (compared 
with  the  KF  performance)  are  collected. 

•  Figures  C.l  through  C.6  present  WF  vs  KF  comparison  for  cases  1  through  6 
respectively.  WF  performance  is  given  for  values  of  w=  0,1,2,3.  The  cases 
correspond  to  following  observer  configurations: 

One  observer 

Two  identical  observers 

Three  identical  observers 

Three   observers,  one  of  whom  3  times  less  accurate  than  the  others 

Five  identical  observers 

Five   observers,  two  of  wnom  3  times  less  accurate  than  the  others 

•  Figures  C.7  and  C.3  present  WF  performance  with  with  different  vaiues  of  w 

Figure  C.7  -  cases  1,2.3 
Figure  C.8  -  cases  4,5,6 

•  Figures  C.9  through  C.14  present  sensitivity  of  WF  performance  with  respect  to 
model  violations.  Figures  C.9  through  C.14  correspond  to  case  1  through  case  6 
respectively. 

•  Figures  C.15  and  C.16  present  comparison  of  WF  and  KF  performances  for 
normally  distributed  measurement  errors.  Comparison  presented  for  cases  la-6a. 
Cases  la-6a  correspond  to  cases  1-6,  but  have  normally  distributed 
measurement  errors  with  the  same  variance  as  Student-t  distributed  errors  in 
original  cases. 

Figure  C.15  presents  cases  la- 3a 

Figure  C.16  presents  cases  4a-6a 

•  Figure  C.17  displays  variability  of  the  simulation  using  different  seeds.  Five 
oairs  of  seeds  are  used.    Variability  is  presented  for  cases  1.2,3.5. 
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Figure  C.2     WF  vs  KF  comparison  for  Case  2  -  Two  Identical  observers. 
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Figure  C.3     WF"  vs  KF  comparison  for  Case  3  -  Three  Identical  observers. 
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Figure  C.4     WF  vs  KF  comparison  for  Case  4  -  Three  Different  observers. 
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Figure  C.5     WF  vs  KF  comparison  for  Case  5  -  Five  Identical  observers. 


44 


CD 

d 


d 


D 


1 

Deviation 

1 
Standard 

. 

— 



KF 
WF 

-u>=0 

— 

WF 

IV —    1 

T=  I 

1 

WF 

u>— 2 

1 

d  = 

i 
i 

1 

3  j  j  3 
1111 

!              ! 

TRAC  LENGTH 

i                !                I 

X 

I 

WF 

1 

i 

I 

1 

0 


20 


40 


60 


SO 


Figure  C.6     WF  vs  KF  comparison  for  Case  0  -  Five  different  observers. 
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Figure  C.9     WF  Sensitivity  Case  1  -  One  observer. 


Gu 


00 


48 


c 


-x- 


WF  correct  assumptions  used 

WF  haif  a"  used 

WF  double  a  used 

WF  haif  r  used 

WF  T— distributed  target  steos 

KF 


"RAC  LENGTH 
i  i  i 


0 


20  40  60  80 

Fieure  CIO    WF  Sensitivity  Case  2  -  Two  Identical  observers. 


100 


■4') 


C\J 


cc 


.VF  correct  assumptions  usee 

-« —  VVF  half  a  used 

-* —  WF  aouDle  d  usee 

- —  WF  half  r  used 

•-X-—  WF  T-distributed  target  steps 
KF 


o 


'RAC  LENGTH 


0 


20 


u 


40 


60 


30 


Figure  C.ll     WF  Sensitivity  Case  3  -  Three  Identical  observers. 
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Figure  C.12     WF  Sensitivity  Case  4  -  Three  DilTcrent  observers. 
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Figure  C.13     WF  Sensitivity  Case  5  •  Five  Identical  observers. 
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Figure  C.14     WF  Sensitivity  Case  6  -  Five  Different  observers. 
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figure  CI  5     WF  vs  KF  comparison  for  Case  la-3a  -  Normal  measurement  errors. 
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Figure  C.16    WF  vs  KF  comparison  for  Case  4a-6a  -  Normal  measurement  errors. 
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